% Estimate effects of US monetary policy shocks using model from
% Antolin-Diaz and Rubio-Ramirez (2018) under different combinations of
% (narrative) sign restrictions. Goal is to: 1) assess sensitivity of
% results to choice of conditional prior; and 2) disentangle
% informativeness of different restrictions.

clear variables
close all
clc

oldFolder = pwd;
cd ..
addpath([oldFolder,'/auxFunctions']);
cd(oldFolder);

%% Import data.
% Data are for Baumeister-Hamilton (2018), 1986.I-2008.III, quarterly
% The variables are:
% 1) Output gap 2) Inflation 3) Federal funds rate (in percent)
load BHdata; data=y; clear y;

% Re-order variables so that the shock of interest is ordered first. 
data = data(:,[end 1:end-1]); % Federal funds rate ordered first

% Create vector of dates in datenum format used for narrative restrictions
dates=zeros(91,1);
dates(12,1)=726438;
dates(33,1)=728326;
dates(20,1)=727168;
dates(52,1)=730029;
dates(62,1)=730942;
dates(68,1)=731521;

% Declare any exogenous variables (other than the constant; leave empty 
% if none).
exog = [];

% Set variable names for plotting figures.
varNames={'Federal funds rate','Output gap','Inflation'};

% Set units for impulse responses. 0 = %, 1 = ppt, 2 = bps.
irfUnits = [2 0 0 0 0 0];

%% Options.
opt.p = 4; % No. of lags in VAR
opt.const = 1; % const = 1 if constant in VAR, = 0 otherwise
opt.jshock = 1; % Index of shock(s) of interest
opt.ivar = [1 2 3];  % Indices of variables of interest
opt.cumIR = []; % Indices of variables for which cumulative impulse response should be used
opt.H = 20; % Maximum horizon for impulse responses
opt.draws = 1000; % No. of draws from posterior of reduced-form parameters (phi) with non-empty identified set
opt.aalpha = 0.68; % Credibility level for credible intervals (0<aalpha<1)
opt.qMax = 100000; % Max no. of draws of orthonormal matrix Q when checking for empty identified set
opt.qDraws = 1000; % No. of draws of Q used to approximate bounds of (conditional) identified set
opt.dispIter = 100; % Print number of draws of phi remaining every dispIter draws
opt.gridLength = 1000; % Size of grid used when computing credible intervals

rng(23032021); % Set seed for random number generator

restr.elastRestr = []; % No elasticity restrictions

%% Sign restrictions only
% Each row of signRestr contains a vector (i,j,h,s,t) representing a
% `traditional' sign restriction, where t is the type of restriction:
% t = 1: the impulse response of the ith variable to the jth shock at the 
% hth horizon is nonnegative (s = 1) or nonpositive (s = -1).
% t = 2: the (ij)th element of A0 is nonnegative (s = 1) or nonpositive 
% (s = -1). 
% Impact sign restrictions:
restr.signRestr = ...
      [1 1 0  1 1; % Response of FFR to monetary policy shock is positive
       2 1 0 -1 1; % Response of output gap to monetary policy shock is negative
       3 1 0 -1 1; % Response of inflation to monetary policy shock is negative
       1 2 0 -1 1; % Response of FFR to cost shock is negative
       2 2 0  1 1; % Response of output gap to cost shock is positive
       3 2 0 -1 1; % Response of inflation to cost shock is negative
       1 3 0  1 1; % Response of FFR to demand shock is positive
       2 3 0  1 1; % Response of output gap to demand shock is positive    
       3 3 0  1 1]; % Response of inflation to demand shock is positive

% restr.signRestr = []; % No sign restrictions on impulse responses

%% Run again with full set of sign and narrative restrictions.
% Each row of shockSignRestr contains a vector (i,s,t) representing the
% shock-sign restriction that the ith shock in period t is >= 0 (s = 1) or
% <= 0 (s = -1). t is equal to the index of the observation.
restr.shockSignRestr = ...
    [1, 1, find(dates == datenum(1988,12,01));
     1, 1, find(dates == datenum(1994,2,01));
     1, -1, find(dates == datenum(1990,12,01));
     1, -1, find(dates == datenum(1998,10,01));
     1, -1, find(dates == datenum(2001,4,01));
     1, -1, find(dates == datenum(2002,11,01))];

%restr.shockSignRestr = []; % No narrative restrictions on sign of shock
restr.hdSignRestr = []; % No restrictions on historical decomposition

% Conduct (robust) Bayesian inference.
fprintf('Estimating BH model with sign and narrative restrictions');

mainfile;

